Romain and I just posted our slides from yesterday's Rcpp workshop and class (preceding the now-ongoing R/Finance conference). You can access the slides via my presentation page, or directly from here as Part 1 (Introduction), Part 2 (Details), Part 3 (Advanced) and Part 4 (Applications). The page also contains a link to the examples as a tar.gz and zip file.

Update: One link corrected.

I realized I never announced this on the blog, so without further ado.... Rcpp Workshop in Chicago on April 28, 2011 This year's R/Finance conference will be preceded by a full-day masterclass on Rcpp and related topics which will be held on Thursday, April 28, 2011, on the University of Illinois at Chicago campus. Join Dirk Eddelbuettel and Romain Francois for six hours of detailed and hands-on instructions and discussions around Rcpp, inline, RInside, RcppArmadillo and other packages---in an intimate small-group setting. The full-day format allows to combine a morning introductory session with a more advanced afternoon session while leaving room for sufficient breaks. There will be about six hours of instructions, a one-hour lunch break and two half-hour coffee breaks. Morning session: "A hands-on introduction to R and C++" The morning session will provide a practical introduction to the Rcpp package (and other related packages). The focus will be on simple and straightforward applications of Rcpp in order to extend R and/or to significantly accelerate the execution of simple functions. The tutorial will cover the inline package which permits embedding of self-contained C, C++ or Fortran code in R scripts. We will also discuss RInside to embed R code in C++ applications, as well as standard Rcpp extension packages such as RcppArmadillo for linear algebra and RcppGSL. Afternoon session: "Advanced R and C++ topics" This afternoon tutorial will provide a hands-on introduction to more advanced Rcpp features. It will cover topics such as writing packages that use Rcpp, how 'Rcpp modules' and the new R ReferenceClasses interact, and how 'Rcpp sugar' lets us write C++ code that is often as expressive as R code. Another possible topic, time permitting, may be writing glue code to extend Rcpp to other C++ projects. We also hope to leave some time to discuss problems brought by the class participants. Prerequisites Knowledge of R as well as general programming knowledge; C or C++ knowledge is helpful but not required. Users should bring a laptop set up so that R packages can be built. That means on Windows, Rtools needs to be present and working, and on OS X the Xcode package should be installed. Registration Registration is available via the R/Finance conference at

http://www.RinFinance.com/register/

or directly at RegOnline

http://www.regonline.com/930153

The cost is USD 500 for the whole day, and space will be limited. Questions Please contact us directly at RomainAndDirk@r-enthusiasts.com.

The RInside package makes it pretty simple and straightforward to embed R, the wonderful statistical programming environment and language, inside of a C++ application. This uses both the robust embedding API provided by R itself, and the higher-level abstractions from our Rcpp package. A number of examples are shown on this blog both here and here as well as on the RInside page; and the source package actually contains well over a dozen complete examples which cover anything from simple examples to parallel use via MPI for parallel computing. Beginning users sometimes ask about how to use RInside inside larger projects. And as I had meant to experiment with embedding inside of the powerful Qt framework anyway, I started to dabble a little. A first result is now in the SVN sources of RInside. My starting point was the classic tkdensity demo that comes with R itself. It is a good point of departure as Tcl/Tk makes it very portable---in fact it should run on every platform that runs R---and quite expressive. And having followed some of the GUI experiments around R over the years, I have also seen various re-implementations using different GUI frameworks. And so I am adding mine to this body of work: The problem I addressed first was actual buildability. For the RInside examples, Romain and I provide a Makefile that just works by making calls to R itself to learn about flags for R, Rcpp and RInside such that all required headers and libraries are found. That is actually relatively straightforward (and documented in our vignettes) but a little intimidating at first---which is why a ready-made Makefile is a good thing. Qt of course uses qmake and the .pro files to encode / resolve dependencies. So task one was to map what our Makefile does into its variables. Turns out that wasn't all that hard:

## -*- mode: Makefile; c-indent-level: 4; c-basic-offset: 4;  tab-width: 8; -*-
##
## Qt usage example for RInside, inspired by the standard 'density
## sliders' example for other GUI toolkits
##
## Copyright (C) 2011  Dirk Eddelbuettel and Romain Francois
TEMPLATE =              app
HEADERS =               qtdensity.h 
SOURCES =               qtdensity.cpp main.cpp
QT +=                   svg
## comment this out if you need a different version of R, 
## and set set R_HOME accordingly as an environment variable
R_HOME =                $$system(R RHOME)
## include headers and libraries for R 
RCPPFLAGS =             $$system($$R_HOME/bin/R CMD config --cppflags)
RLDFLAGS =              $$system($$R_HOME/bin/R CMD config --ldflags)
RBLAS =                 $$system($$R_HOME/bin/R CMD config BLAS_LIBS)
RLAPACK =               $$system($$R_HOME/bin/R CMD config LAPACK_LIBS)
## if you need to set an rpath to R itself, also uncomment
#RRPATH =               -Wl,-rpath,$$R_HOME/lib
## include headers and libraries for Rcpp interface classes
RCPPINCL =              $$system($$R_HOME/bin/Rscript -e \'Rcpp:::CxxFlags\(\)\')
RCPPLIBS =              $$system($$R_HOME/bin/Rscript -e \'Rcpp:::LdFlags\(\)\')
## for some reason when building with Qt we get this each time
## so we turn unused parameter warnings off
RCPPWARNING =           -Wno-unused-parameter 
## include headers and libraries for RInside embedding classes
RINSIDEINCL =           $$system($$R_HOME/bin/Rscript -e \'RInside:::CxxFlags\(\)\')
RINSIDELIBS =           $$system($$R_HOME/bin/Rscript -e \'RInside:::LdFlags\(\)\')
## compiler etc settings used in default make rules
QMAKE_CXXFLAGS +=       $$RCPPWARNING $$RCPPFLAGS $$RCPPINCL $$RINSIDEINCL
QMAKE_LFLAGS +=         $$RLDFLAGS $$RBLAS $$RLAPACK $$RCPPLIBS $$RINSIDELIBS
## addition clean targets
QMAKE_CLEAN +=          qtdensity Makefile

The double dollar signs and escaping of parentheses are a little tedious, but hey it works and expands the compiler and linker flags such that everything <emjust works="Works">. The code itself is pretty straightforward too. We instantiate the RInside object as well as the main Qt application object. We then instantiate a new object of class QtDensity that will launch the main widget; it is given a reference to the RInside object.

// -*- mode: C++; c-indent-level: 4; c-basic-offset: 4;  tab-width: 8; -*-
//
// Qt usage example for RInside, inspired by the standard 'density
// sliders' example for other GUI toolkits
//
// Copyright (C) 2011  Dirk Eddelbuettel and Romain Francois
#include <QApplication>
#include "qtdensity.h"
int main(int argc, char *argv[])
 
    RInside R(argc, argv);  		// create an embedded R instance
    QApplication app(argc, argv);
    QtDensity qtdensity(R);
    return app.exec();
 

The definition of the main object is pretty simple: a few private variables, and a few functions to interact with the GUI and get values from the radio buttons, slider or input field---as well as functions to update the chart or re-draw the random variables.

// -*- mode: C++; c-indent-level: 4; c-basic-offset: 4;  tab-width: 8; -*-
//
// Qt usage example for RInside, inspired by the standard 'density
// sliders' example for other GUI toolkits
//
// Copyright (C) 2011  Dirk Eddelbuettel and Romain Francois
#ifndef QTDENSITY_H
#define QTDENSITY_H
#include <RInside.h>
#include <QMainWindow>
#include <QHBoxLayout>
#include <QSlider>
#include <QSpinBox>
#include <QLabel>
#include <QTemporaryFile>
#include <QSvgWidget>
class QtDensity : public QMainWindow
 
    Q_OBJECT
public:
    QtDensity(RInside & R);
private slots:
    void getBandwidth(int bw);
    void getKernel(int kernel);
    void getRandomDataCmd(QString txt);
    void runRandomDataCmd(void);
private:
    void setupDisplay(void);    // standard GUI boilderplate of arranging things
    void plot(void);            // run a density plot in R and update the
    void filterFile(void);      // modify the richer SVG produced by R
    QSvgWidget *m_svg;          // the SVG device
    RInside & m_R;              // reference to the R instance passed to constructor
    QString m_tempfile;         // name of file used by R for plots
    QString m_svgfile;          // another temp file, this time from Qt
    int m_bw, m_kernel;         // parameters used to estimate the density
    QString m_cmd;              // random draw command string
 ;
#endif

Lastly, no big magic in the code either (apart from the standard magic provided by RInside). A bit of standard GUI layouting, and then some functions to pick values from the inputs as well as to compute / update the output. One issue is worth mentioning. The screenshot and code show the second version of this little application. I built a first one using a standard portable network graphics (png) file. That was fine, but not crisp as png is a pixel format so I went back and experimented with scalable vector graphics (svg) instead. One can create svg output with R in a number of ways, one of which is the cairoDevice package by Michael Lawrence (who also wrote RGtk2 and good chunks of Ggobi). Now, it turns out that Qt displays the so-called SVG tiny standard whereas R creates a fuller SVG format. Some discussion with Michael reveals that one can modify the svg file suitably (which is what the function filterFile below does) and it all works. Well: almost. There is a bug (and Michael thinks it is the SVG rendering) in which the density estimate does not get clipped to the plotting region.

// -*- mode: C++; c-indent-level: 4; c-basic-offset: 4;  tab-width: 8; -*-
//
// Qt usage example for RInside, inspired by the standard 'density
// sliders' example for other GUI toolkits -- this time with SVG
//
// Copyright (C) 2011  Dirk Eddelbuettel and Romain Francois
#include <QtGui>
#include "qtdensity.h"
QtDensity::QtDensity(RInside & R) : m_R(R)
 
    m_bw = 100;                 // initial bandwidth, will be scaled by 100 so 1.0
    m_kernel = 0;               // initial kernel: gaussian
    m_cmd = "c(rnorm(100,0,1), rnorm(50,5,1))"; // simple mixture
    m_R["bw"] = m_bw;           // pass bandwidth to R, and have R compute a temp.file name
    m_tempfile = QString::fromStdString(Rcpp::as<std::string>(m_R.parseEval("tfile <- tempfile()")));
    m_svgfile = QString::fromStdString(Rcpp::as<std::string>(m_R.parseEval("sfile <- tempfile()")));
    m_R.parseEvalQ("library(cairoDevice)");
    setupDisplay();
 
void QtDensity::setupDisplay(void)   
    QWidget *window = new QWidget;
    window->setWindowTitle("Qt and RInside demo: density estimation");
    QSpinBox *spinBox = new QSpinBox;
    QSlider *slider = new QSlider(Qt::Horizontal);
    spinBox->setRange(5, 200);
    slider->setRange(5, 200);
    QObject::connect(spinBox, SIGNAL(valueChanged(int)), slider, SLOT(setValue(int)));
    QObject::connect(slider, SIGNAL(valueChanged(int)), spinBox, SLOT(setValue(int)));
    spinBox->setValue(m_bw);
    QObject::connect(spinBox, SIGNAL(valueChanged(int)), this, SLOT(getBandwidth(int)));
    QLabel *cmdLabel = new QLabel("R command for random data creation");
    QLineEdit *cmdEntry = new QLineEdit(m_cmd);
    QObject::connect(cmdEntry,  SIGNAL(textEdited(QString)), this, SLOT(getRandomDataCmd(QString)));
    QObject::connect(cmdEntry,  SIGNAL(editingFinished()), this, SLOT(runRandomDataCmd()));
    QGroupBox *kernelRadioBox = new QGroupBox("Density Estimation kernel");
    QRadioButton *radio1 = new QRadioButton("&Gaussian");
    QRadioButton *radio2 = new QRadioButton("&Epanechnikov");
    QRadioButton *radio3 = new QRadioButton("&Rectangular");
    QRadioButton *radio4 = new QRadioButton("&Triangular");
    QRadioButton *radio5 = new QRadioButton("&Cosine");
    radio1->setChecked(true);
    QVBoxLayout *vbox = new QVBoxLayout;
    vbox->addWidget(radio1);
    vbox->addWidget(radio2);
    vbox->addWidget(radio3);
    vbox->addWidget(radio4);
    vbox->addWidget(radio5);
    kernelRadioBox->setMinimumSize(260,140);
    kernelRadioBox->setMaximumSize(260,140);
    kernelRadioBox->setSizePolicy(QSizePolicy::Fixed, QSizePolicy::Fixed);
    kernelRadioBox->setLayout(vbox);
    QButtonGroup *kernelGroup = new QButtonGroup;
    kernelGroup->addButton(radio1, 0);
    kernelGroup->addButton(radio2, 1);
    kernelGroup->addButton(radio3, 2);
    kernelGroup->addButton(radio4, 3);
    kernelGroup->addButton(radio5, 4);
    QObject::connect(kernelGroup, SIGNAL(buttonClicked(int)), this, SLOT(getKernel(int)));
    m_svg = new QSvgWidget();
    runRandomDataCmd();         // also calls plot()
    QGroupBox *estimationBox = new QGroupBox("Density estimation bandwidth (scaled by 100)");
    QHBoxLayout *spinners = new QHBoxLayout;
    spinners->addWidget(spinBox);
    spinners->addWidget(slider);
    QVBoxLayout *topright = new QVBoxLayout;
    topright->addLayout(spinners);
    topright->addWidget(cmdLabel);
    topright->addWidget(cmdEntry);
    estimationBox->setMinimumSize(360,140);
    estimationBox->setMaximumSize(360,140);
    estimationBox->setSizePolicy(QSizePolicy::Fixed, QSizePolicy::Fixed);
    estimationBox->setLayout(topright);
    QHBoxLayout *upperlayout = new QHBoxLayout;
    upperlayout->addWidget(kernelRadioBox);
    upperlayout->addWidget(estimationBox);
    QHBoxLayout *svglayout = new QHBoxLayout;
    svglayout->addWidget(m_svg);
    QVBoxLayout *outer = new QVBoxLayout;
    outer->addLayout(upperlayout);
    outer->addLayout(svglayout);
    window->setLayout(outer);
    window->show();
 
void QtDensity::plot(void)  
    const char *kernelstrings[] =   "gaussian", "epanechnikov", "rectangular", "triangular", "cosine"  ;
    m_R["bw"] = m_bw;
    m_R["kernel"] = kernelstrings[m_kernel]; // that passes the string to R
    std::string cmd1 = "Cairo(width=6,height=6,pointsize=10,surface='svg',filename=tfile); "
                       "plot(density(y, bw=bw/100, kernel=kernel), xlim=range(y)+c(-2,2), main=\"Kernel: ";
    std::string cmd2 = "\"); points(y, rep(0, length(y)), pch=16, col=rgb(0,0,0,1/4));  dev.off()";
    std::string cmd = cmd1 + kernelstrings[m_kernel] + cmd2; // stick the selected kernel in the middle
    m_R.parseEvalQ(cmd);
    filterFile();               // we need to simplify the svg file for display by Qt 
    m_svg->load(m_svgfile);
 
void QtDensity::getBandwidth(int bw)  
    if (bw != m_bw)  
        m_bw = bw;
        plot();
     
 
void QtDensity::getKernel(int kernel)  
    if (kernel != m_kernel)  
        m_kernel = kernel;
        plot();
     
 
void QtDensity::getRandomDataCmd(QString txt)  
    m_cmd = txt;
 
void QtDensity::runRandomDataCmd(void)  
    std::string cmd = "y <- " + m_cmd.toStdString();
    m_R.parseEvalQ(cmd);
    plot();                     // after each random draw, update plot with estimate
 
void QtDensity::filterFile()  
    // cairoDevice creates richer SVG than Qt can display
    // but per Michaele Lawrence, a simple trick is to s/symbol/g/ which we do here
    QFile infile(m_tempfile);
    infile.open(QFile::ReadOnly);
    QFile outfile(m_svgfile);
    outfile.open(QFile::WriteOnly   QFile::Truncate);
    
    QTextStream in(&infile);
    QTextStream out(&outfile);
    QRegExp rx1("<symbol"); 
    QRegExp rx2("</symbol");    
    while (!in.atEnd())  
        QString line = in.readLine();
        line.replace(rx1, "<g"); // so '<symbol' becomes '<g ...'
        line.replace(rx2, "</g");// and '</symbol becomes '</g'
        out << line << "\n";
     
    infile.close();
    outfile.close();
 

What the little application does is actually somewhat neat for the few lines. One key features is that the generated data can be specified directly by an R expression which allows for mixtures (as shown, and as is the default). With that it easy to see how many points are needed in the second hump to make the estimate multi-modal, and how much of a distance between both centers is needed and so on. Obviously, the effect of the chosen kernel and bandwidth can also be visualized. And with the chart the being a support vector graphics display, we can resize and scale at will and it still looks crisp. The code (for both the simpler png variant and the svg version shown here) is in the SVN repository for RInside and will be in the next release. Special thanks to Michael Lawrence for patiently working through some svg woes with me over a few emails.

Update: Some typos fixed.

Update 2: Two URLs corrected.

My friend Jeff deserves a sincere congratulation for finally unveiling his rebranded R consultancy Lemnica. One notable feature of the new website is a section called esoteric R which discusses less frequently-visited corners of the R world. It even boasts its own CRAN package esotericR with the example sources. esoteric R currently holds two articles. Jeff had sent me the one about introducing closures a while back, and I like it and may comment at another time. What caught me by surprise when Lemnica finally opened was the other article: R calling C. It is a fine article motivated by all the usual reasons that are e.g. mentioned in the Google Tech Talk which Romain and I gave last October about our work around Rcpp. But it is just not simple. Allow me to explain. When Jeff showed this C language file

#include <R.h>
#include <Rinternals.h>
SEXP esoteric_rev (SEXP x)  
  SEXP res;
  int i, r, P=0;
  PROTECT(res = allocVector(REALSXP, length(x))); P++;
  for(i=length(x), r=0; i>0; i--, r++)  
     REAL(res)[r] = REAL(x)[i-1];
   
  copyMostAttrib(x, res);
  UNPROTECT(P);
  return res;
 

and then needs several paragraphs to explain what is going on, what is needed to compile and then how to load it --- I simply could not resist. Almost immediately, I emailed back to him something as simple as this using both our Rcpp package as well as the wonderful inline package by Oleg which Romain and I more or less adopted:

library(inline)  ## for cxxfunction()
src <- 'Rcpp::NumericVector x = Rcpp::NumericVector(xs);
        std::reverse(x.begin(), x.end());
        return(x);'
fun <- cxxfunction(signature(xs="numeric"), body=src, plugin="Rcpp")
fun( seq(0, 1, 0.1) )

Here we load inline, and then define a three-line C++ program using facilities from our Rcpp package. All we need to revert a vector is to first access its R object in C++ by instantiating the R vector as a NumericVector. These C++ classes then provide iterators which are compatible with the Standard Template Library (STL). So we simply call the STL function reverse pointing the beginning and end of the vector, and are done! Rcpp then allows us the return the C++ vector which it turns into an R vector. Efficient in-place reversal, just like Jeff had motivated, in three lines. Best of all, we can execute this from within R itself:

R> library(inline)  ## for cxxfunction()
R> src <- 'Rcpp::NumericVector x = Rcpp::NumericVector(xs);
+         std::reverse(x.begin(), x.end());
+         return(x);'
R> fun <- cxxfunction(signature(xs="numeric"), body=src, plugin="Rcpp")
R> fun( seq(0, 1, 0.1) )
 [1] 1.0 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0.0
R>

Lastly, Jeff shows a more complete example wherein a new vector is created, and any potential attributes are copied as well. Naturally, we can do that too. First, we used clone() to make a deep copy (ie forcing creation of a new object rather than a mere proxy) and use the same R API function he accessed---but it our case both prefixed with ::Rf_ for R remapping (to protect clashed with other functions with identical names) and a global namespace identifier (as it is a global C function from R).

R> library(inline)
R> src <- 'Rcpp::NumericVector x = Rcpp::clone<Rcpp::NumericVector>(xs);
+         std::reverse(x.begin(), x.end());
+         ::Rf_copyMostAttrib(xs, x);
+         return(x);'
R> fun <- cxxfunction(signature(xs="numeric"), body=src, plugin="Rcpp")
R> obj <- structure(seq(0, 1, 0.1), obligatory="hello, world!")
R> fun(obj)
 [1] 1.0 0.9 0.8 0.7 0.6 0.5 0.4 0.3 0.2 0.1 0.0
attr(,"obligatory")
[1] "hello, world!"
R> obj
 [1] 0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1.0
attr(,"obligatory")
[1] "hello, world!"
R>

Both the obj variable and the new copy contain the desired data attribute, the new copy is reversed, the original is untouched---and all in four lines of C++ called via one inline call. I have now been going on for over one hundred lines yet I never had to mention memory management, pointers, PROTECT or other components of the R API for C. Hopefully, this short writeup provided an idea of why Romain and I think Rcpp is the way to go for creating C/C++ functions for extending and enhancing R.

I've decided to take over the maintenance of the unofficial emacs-snapshot Debian packages that were maintained by Romain Francoise. They are available on a dedicated page.

The text below went out as a post to the r-packages list a few days ago, but I thought it would make sense to post it on the blog too. So with a little html markup... Summary Version 0.9.0 of the Rcpp package is now on CRAN and its mirrors. This release marks another step in the development of the package, and a few key points are highlighted below. More details are in the NEWS and ChangeLog files included in the package. Overview Rcpp is an R package and associated C++ library that facilitates integration of C++ code in R packages. The package features a complete set of C++ classes (Rcpp::IntegerVector, Rcpp:NumericVector, Rcpp::Function, Rcpp::Environment, ...) that makes it easier to manipulate R objects of matching types (integer vectors, functions, environments, etc ...). Rcpp takes advantage of C++ language features such as the explicit constructor / destructor lifecycle of objects to manage garbage collection automatically and transparently. We believe this is a major improvement over use of PROTECT / UNPROTECT. When an Rcpp object is created, it protects the underlying SEXP so that the garbage collector does not attempt to reclaim the memory. This protection is withdrawn when the object goes out of scope. Moreover, users generally do not need to manage memory directly (via calls to new / delete or malloc / free) as this is done by the Rcpp classes or the corresponding STL containers. A few key points about Rcpp:

• a rich API covering all core R data types including vectors, matrices, functions, environments, ... (with the exeception of factors which are less useful in C++)
• seamless (bi-directional) data interchange between R and C++
• possibility of inline use permitting definition, compilation, linking and loading of C++ functions directly from R
• extensive documentation now covering eight vignettes
• exception handling and error propagation back to R
• extensive test suite using RUnit covering over 700 tests
• extension packages RcppArmadillo and RcppGSL provide easy-to-use integration with the Armadillo (linear algebra) and GNU GSL librasries
• increasing adoption among R users and package developers with now twenty packages from CRAN or BioConductor depending on Rcpp
• support for the legacy 'classic' Rcpp is now provided by the RcppClassic package which is being released concurrently with Rcpp 0.9.0

Several key features were added during the 0.8.* cycles and are described below. Rcpp sugar Rcpp now provides syntactic sugar: vectorised expressions at the C++ level which are motivated by the corresponding R expressions. This covers operators (binary arithmetic, binary logical, unary), functions (producing single logical results, mathematical functions and d/p/q/r statistical functions). Examples comprises anything from ifelse() to pmin()/pmax() or A really simply example is a function

    SEXP foo( SEXP xx, SEXP yy) 
        NumericVector x(xx), y(yy) ;
        return ifelse( x < y, x*x, -(y*y) ) ;
     

which deploys the sugar 'ifelse' function modeled after the corresponding R function. Another simple example is

    double square( double x) 
        return x*x ;
     
    SEXP foo( SEXP xx ) 
        NumericVector x(xx) ;
        return sapply( x, square ) ;
     

where use the sugar function 'sapply' to sweep a simple C++ function which operates elementwise across the supplied vector. The Rcpp-sugar vignette describes sugar in more detail. Rcpp modules Rcpp modules are inspired by Boost.Python and make exposing C++ functions or classes to R even easier. A first illustration is provided by this simple C++ code snippet

    const char* hello( const std::string& who ) 
        std::string result( "hello " ) ;
        result += who ;
        return result.c_str() ;
     
    RCPP_MODULE(yada) 
        using namespace Rcpp ;
        function( "hello", &hello ) ;
     

which (after compiling and loading) we can access in R as

    yada <- Module( "yada" )
    yada$hello( "world" )

In a similar way, C++ classes can be exposed very easily. Rcpp modules are also described in more detail in their own vignette. Reference Classes R release 2.12.0 introduced Reference Classes. These are formal S4 classes with the corresponding dispatch method, but passed by reference and easy to use. Reference Classes can also be exposed to R by using Rcpp modules. Extension packackages The RcppArmadillo package permits use of the advanced C++ library 'Armadillo, a C++ linear algebra library aiming towards a good balance between speed and ease of use, providing integer, floating point and complex matrices and vectors with lapack / blas support via R. Armadillo uses templates for a delayed evaluation approach is employed (during compile time) to combine several operations into one and reduce (or eliminate) the need for temporaries. Armadillo is useful if C++ has been decided as the language of choice, rather than another language like Matlab or Octave, and aims to be as expressive as the former. Via Rcpp and RcppArmadillo, R users now have easy access to this functionality. Examples are provided in the RcppArmadillo package. The RcppGSL package permits easy use of the GNU Scientific Library (GSL), a collection of numerical routines for scientifc computing. It is particularly useful for C and C++ programs as it provides a standard C interface to a wide range of mathematical routines such as special functions, permutations, combinations, fast fourier transforms, eigensystems, random numbers, quadrature, random distributions, quasi-random sequences, Monte Carlo integration, N-tuples, differential equations, simulated annealing, numerical differentiation, interpolation, series acceleration, Chebyshev approximations, root-finding, discrete Hankel transforms physical constants, basis splines and wavelets. There are over 1000 functions in total with an extensive test suite. The RcppGSL package provides an easy-to-use interface between GSL data structures and R using concepts from Rcpp. The RcppGSL package also contains a vignette with more documentation. Legacy 'classic' API Packages still using code interfacing the initial 'classic' Rcpp API are encouraged to migrate to the new API. Should a code transition not be possible, backwards compatibility is provided by the RcppClassic package released alongside Rcpp 0.9.0. By including RcppClassic.h and building against the RcppClassic package and library, vintage code can remain operational using the classic API. The short vignette in the RcppClassic package has more details. Documentation The package contains a total of eight vignettes the first of which provides a short and succinct introduction to the Rcpp package along with several motivating examples. Links

• Rcpp main page
• R-forge project page
• Dirk's blog section about Rcpp
• Romain's blog section about Rcpp

Support Questions about Rcpp should be directed to the Rcpp-devel mailing list https://lists.r-forge.r-project.org/cgi-bin/mailman/listinfo/rcpp-devel

Dirk Eddelbuettel, Romain Francois, Doug Bates and John Chambers
December 2010

A new package RcppDE has been uploaded in a first version 0.1.0 to CRAN. It provides differential evolution optimisation---a variant of stochastic optimisation that is similar to genetic algorithms but particularly suitable for the floating-point representations common in numerical optimisation. It builds of on the nice DEoptim package by Ardia et al, but reimplements the algorithm in C++ (rather than C) using a large serving of Rcpp and RcppArmadillo. I worked on this on for a few evenings and weekends in October and November and then spent a few more evenings writing a paper / vignette (which is finished as a very first draft now) about it. This was an interesting and captivating problem as I had worked on genetic algorithms going back quite some time to the beginning and then again the end of graduate school (and traces of that early work are near the bottom of my presentations page). So what got me started? DEoptim is a really nice package, but it is implemented in old-school C. There is nothing wrong with that per se, but at the same time that I was wrestling with GAs, I also taught myself C++ which, to put it simply, offers a few more choices to the programmer. I like having those choices. And with all the work that Romain and I have put into Rcpp, I was curious how far I could push this cart if I were to move it along. I made a bet with myself starting from the old saw shorter, easier, faster: pick any two. Would it be possible to achieve all three of these goals? DEoptim, and I take version 2.0-7 as my reference point here, is pretty efficiently yet verbosely coded. Copying a vector takes a loop with an assignment for each element, copying a matrix does the same using two loops. Replacing that with a single statement in C++ is pretty easy. We also have a few little optimisations behind the scenes here and there in Rcpp: would all that be enough to move the needle in terms of performance? And the same time, DEoptim is also full of the uses of the old R API which we often point to in the Rcpp documentation so fixing readibility should be a relatively low-hanging fruit. To cut a long story short, I was able to reduce code size quite easily by using a combination of C++ and Rcpp idioms. I was also able to get to faster: the paper / vignette demostrates consistent speed improvements on all setups that I tested (three standard functions on three small and three larger parameter vectors). More important speed gains were achieved by allowing use of objective functions that are written in C++ which again is both possible and easy thanks to Rcpp. That leaves easier to prove: adding compiled objective functions is one indication; further proof could be provided by, say, moving the inner loop to parallel execution thanks to Open MP which I may attempt over the next few months. So far I'd like to give myself about half a point here. So not quite yet shorter, easier, faster: pick any three, but working on it. Over the next few days I may try to follow up with a blog post or two contrasting some code examples and maybe showing a chart from the vignette.

Romain pushed verion 0.3.8 of inline to CRAN earlier today, and I just updated the Debian package. This version adds an internal performance enhancement which is obtained by making due with fewer reads. The short NEWS file entry follows:

0.3.8   2010-12-07
    o   faster cfunction and cxxfunction by loading and resolving the routine
        at "compile" time

Earlier in the year, Romain and I did a bunch of initial work on a wrapper from R to the GNU GSL by way of our Rcpp package for seamless R and C++ integration. But other work kept us busy and this fell a little to the side. We have now found some time to finish this work for a first release, together with a nicely detailed eleven page package vignette. As of today, the package is now a CRAN package, and Romain already posted a nice announcement on his blog and on the rcpp-devel list. So what does RcppGSL do? I gave the package its own webpage here as well and listed these points as key features of RcppGSL:

• templated vector and matrix classes: these are similar to Rcpp's own vector and matrix classes, but really are just smart pointers around the C structure expected by the library
• this means you can transfer data from R to your GSL-using programs in pretty much the same way you would in other C++ programs using Rcpp---by relying on the Rcpp::as() and Rcpp::wrap() converterrs
• at the C++ level you can use these GSL vectors in a more C++-alike way (using eg foo[i] to access an element at index i)
• yet at the same time you can just pass these vector and matrix objects to the GSL functions expecting its C objects: thanks to some cleverness in these classes they pass the right object on (see the example below)
• we also provide the lightweight views for vectors and matrices as the GSL API uses these in many places.

Also provided is a simple example which is a simple implementation of a column norm (which we could easily compute directly in R, but we are simply re-using an example from Section 8.4.14 of the GSL manual):

#include <RcppGSL.h>
#include <gsl/gsl_matrix.h>
#include <gsl/gsl_blas.h>
extern "C" SEXP colNorm(SEXP sM)  
  try  
        RcppGSL::matrix<double> M = sM;     // create gsl data structures from SEXP
        int k = M.ncol();
        Rcpp::NumericVector n(k);           // to store results
        for (int j = 0; j < k; j++)  
            RcppGSL::vector_view<double> colview = gsl_matrix_column (M, j);
            n[j] = gsl_blas_dnrm2(colview);
         
        M.free() ;
        return n;                           // return vector
    catch( std::exception &ex )  
        forward_exception_to_r( ex );
    catch(...)  
        ::Rf_error( "c++ exception (unknown reason)" );
   
  return R_NilValue; // -Wall
 

This example function is implemented in an example package contained in the RcppGSL package itself -- so that users have a complete stanza to use in their packages. This will then build a user package on Linux, OS X and Windows provided the GSL is installed (and on Windows you have to do all the extra steps of defining an environment variable pointing to and of course install Rtools to build in the first place---Linux and OS X are so much easier for development). Another complete example is in the package itself and provides a faster (compiled) alternative to the standard lm() function in R; this example is the continuation of the same example I had in several versions of my Intro to HPC with R tutorials and in the Rcpp package itself as an early example. We will try to touch base with CRAN package authors using both GSL and Rcpp to see how this can help them. The API in our package may well be incomplete, but we are always happy to try to respond to requests for additional features brought to our attention, preferably via the rcpp-devel list. More information is on the RcppGSL page. Questions, comments etc should go to the rcpp-devel mailing list off the R-Forge page.

A bug-fix release 0.8.8 of Rcpp is now available. It is awaiting processing at CRAN, and will be uploaded to Debian once processed at CRAN. In the meantime, sources are available from my local directory here. This release follows on the heels of 0.8.7, but contains fixes for a few small things Romain and I had noticed over the last two weeks since releasing 0.8.7 and contains only a small number of new tweaks. The NEWS entry follows below:

0.8.8   2010-11-01
    o   New syntactic shortcut to extract rows and columns of a Matrix. 
        x(i,_) extracts the i-th row and x(_,i) extracts the i-th column. 
    
    o   Matrix indexing is more efficient. However, faster indexing is
        disabled if g++ 4.5.0 or later is used.
    o   A few new Rcpp operators such as cumsum, operator=(sugar)
    o   Variety of bug fixes:
        - column indexing was incorrect in some cases
        - compilation using clang/llvm (thanks to Karl Millar for the patch)
        - instantation order of Module corrected
        - POSIXct, POSIXt now correctly ordered for R 2.12.0 

As always, even fuller details are on the Rcpp Changelog page and the Rcpp page which also leads to the downloads, the browseable doxygen docs and zip files of doxygen output for the standard formats. A local directory has source and documentation too. Questions, comments etc should go to the rcpp-devel mailing list off the R-Forge page

Last Friday, Romain and I were guests of the R intergrouplet (what an adorable name!) at Google's headquarter in Mountain View. This arose out of discussions following useR! 2010 where we met Google's Murray Stokely. There appears to be ever increasing use of R at Google, and so it was a great opportunity to give a Google Tech Talk about R and C++ integration --- centered around our Rcpp, RInside and RProtoBuf packages which facilitate interoperability between R and C++. A video recording of our ninety-minute talk is already available via the YouTube channel for Google Tech Talks. The (large) pdf with slides (which Romain had already posted on slideshare) is also available from my presentations page. The remainder of the weekend was nice too (with the notably exception of the extremly sucky weather). We got to to spend some time at the Google Summer of Code Mentor Summit which is always a fun event and a great way to meet other open source folks in person. And we also took one afternoon off to spend some with John Chambers discussing further work involving Rcpp and the new ReferenceClasses that appeared in the just-released R version 2.12.0. This should be a nice avenue to further integrate R and C++ in the near future.

With the scheduled R release of version 2.12.0 this morning, we have just uploaded version 0.8.7 of Rcpp to CRAN; Debian will follow shortly once the autobuilders have processed R 2.12.0 This Rcpp release depends on R 2.12.0 as two things have changed. First, we play along with change in R concerning the ordering of inheritance for time classes. But secondly, and more importantly, we support in Rcpp the corresponding change R itself which brings the new ReferenceClasses. Here is corresponding bit from R's NEWS file for R 2.12.0:

    o A facility for defining reference-based S4 classes (in the OOP
      style of Java, C++, etc.) has been added experimentally to
      package methods; see ?ReferenceClasses.
[...]
    o An experimental new programming model has been added to package
      methods for reference (OOP-style) classes and methods.  See
      ?ReferenceClasses.

This was made possible in large part by code committed by John Chambers (whom we had welcomed recently as a co-author to Rcpp) building on the changes he made to R 2.12.0 itself, as well on the work Romain had done with 'Rcpp Modules'. The R help page for ReferenceClasses carries a reference (bad pun) to Rcpp 0.8.7 so these two releases do go together. This should be a lot of fun over the next little while: S3, S4, and now ReferenceClasses. We also made a number of internal changes some of which leads to speed-ups and internal improvement. The NEWS entry follows below:

0.8.7   2010-10-15
    o   As of this version, Rcpp depends on R 2.12 or greater as it interfaces 
        the new reference classes (see below) and also reflects the POSIXt class
        reordering both of which appeared with R version 2.12.0
    o   new Rcpp::Reference class, that allows internal manipulation of R 2.12.0
        reference classes. The class exposes a constructor that takes the name
        of the target reference class and a field(string) method that implements
        the proxy pattern to get/set reference fields using callbacks to the 
        R operators "$" and "$<-" in order to preserve the R-level encapsulation
    o   the R side of the preceding item allows methods to be written
        in R as per ?ReferenceClasses, accessing fields by name and
        assigning them using "<<-".  Classes extracted from modules
        are R reference classes.  They can be subclassed in R, and/or R methods
        can be defined using the $methods(...) mechanism.
    o   internal performance improvements for Rcpp sugar as well as an added
        'noNA()' wrapper to omit tests for NA values -- see the included
        examples in inst/examples/convolveBenchmarks for the speedups
    o   more internal performance gains with Functions and Environments

As always, even fuller details are in Rcpp Changelog page and the Rcpp page which also leads to the downloads, the browseable doxygen docs and zip files of doxygen output for the standard formats. A local directory has source and documentation too. Questions, comments etc should go to the rcpp-devel mailing list off the R-Forge page

---

A brand new and shiny release of RProtoBuf, now at version 0.2.0, arrived on CRAN earlier today. RProtoBuf provides GNU R bindings for the Google Protobuf data encoding library used and release by Google and others. This is only the second release after 0.1-0 more than six months ago. Given that Rcpp is such a key ingedrient for RProtoBuf, and that Rcpp underwent so many exciting changes itself, Romain and I never got around to releasing new versions of RProtoBuf. This version is now much closer to the actual C++ API and fairly feature rich. We summarised a few of these new things in the presentation at useR! 2010. There is more information at the RProtoBuf page; there is a draft package vignette, a 'quick' overview vignette and a unit test summary vignette. Questions, comments etc should go to the rprotobuf mailing list off the RProtoBuf page at R-Forge.

---

Back in August 2007 I looked at the state of embedded md5sums in Debian packages and found that approximately 3% of the files in the archive didn't have checksums. Three years later, things have improved: only 0.76% of the archive is now missing checksums (sid, main/contrib/non-free). (See this lintian report for the list of affected packages.) Since then there's also been various discussions on this subject and there is now a policy bug open to make md5sums a requirement (at the "should" level). There is also a wishlist bug against debhelper to turn dh_md5sums into dh_checksums with a stronger hash algorithm, but MD5 still being good enough for simple integrity checking, it seems rather pointless to upgrade the algorithm without a trust path in the form of in-package signatures ala RPM... Anyway, what's perhaps more surprising is the growth of the distribution in only three years: sid has gone from 20774 to 30314 packages, a 45% increase. Similarly, the number of regular files has gone from approximately 2 million to just above 2.9 million. Indeed, looking at our last five releases, the distribution's growth is impressive:

• woody (2002): 8273 packages
• sarge (2005): 15195 packages (+83.7%)
• etch (2007): 18043 packages (+18.7%)
• lenny (2009): 22277 packages (+23.5%)
• squeeze (2010?): 28870 packages (+29.6%)

Whether or not that is a good thing is, of course, yet to be determined. As a data point, I used the UDD to know how many of these thousands of packages are actually used, and to my surprise, 22321 binary packages have a popcon installation count that is less than 500! (By comparison, dpkg's installation count is 89393.) So while each new release adds lots of packages, the majority of them have very few users. (If you want to check yourself, the query I used is select p.package, version, insts from packages p, popcon where (p.architecture = 'i386' or p.architecture = 'all') and p.release = 'squeeze' and p.package = popcon.package and popcon.insts < 500 order by insts;.)